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Abstract 

A zonal approach for direct computation of sound 
generation and propagation from a supersonic jet is in- 
vestigated. The present work splits the computational 
domain into a non-linear, acoustic-source regime and a 
linear acoustic wave propagation regime. In the non-lin- 
ear regime, the unsteady flow is governed by the large- 
scale equations, which are the filtered compressible 
Navier-Stokes equations. In the linear acoustic regime, 
the sound wave propagation is described by the linear- 
ized Euler equations. Computational results are present- 
ed for a supersonic jet at M=2.1. It is demonstrated that 
no spurious modes are generated in die matching region 
and the computational expense is reduced substantially 
as opposed to fully large-scale simulation. 

1. Introduction 

In theory, direct numerical simulation (DNS) based 
on the compressible Navier-Stokes equations provide 
both the flow fluctuadons and the acoustic field Howev- 
er, die resolution requirement for high-Reynolds number 
turbulent flows makes direct numerical simulation im- 
practical due to current computer limitations. It is known 
that the large scale structure is responsible for the gener- 
ation of the dominant part of supersonic jet noise [1-6]. 
This indicates that it is appropriate to perform large-eddy 
simulations(LES) to accurately capture the large scales 
of motion while modelling the sub-grid scale turbulence. 

The use of large-eddy simulations (LES) as a tool 
for prediction of the jet noise source has been proposed 
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by Mankbadi et aL [7], in which, the computed fluctuat- 
ing sound source in the near field was then used to obtain 
the far field sound through the application of LighthilTs 
theory. Due to the non-compactness of the source and the 
need to account for the retarded time, accurate applica- 
tion of LighthilTs theory requires prohibitive computer 
storage. Later, Mankbadi et al. [8] extended the compu- 
tational domain to the acoustic field, where large-scale 
equations are solved for both the sound source and the 
acoustic fields for the axisymmetric case. The approach 
was extended in Shih et al. [9] to the three dimensional 
case, but was restricted to the near region due to comput- 
er limitations. An alternative approach for extending the 
three dimensional computation to the far-field is needed. 
This issue is addressed in this paper. 

An acoustic calculation can be viewed as consisting 
of two parts, one describing the non-linear generation of 
sound, the other describing the linear propagation of 
sound. All non-linear flow effects and source generation 
are confined to the near field, and can be computed by 
large-scale equations. There are several approaches to 
calculate the sound propagation once the source has been 
identified, such as acoustic analogy, Kirchhoff s method 
[10], and linearized Euler equations. Freund et al. [11] 
have studied the matching of near/far-field equation sets 
for computations of aerodynamic sound. Their results of 
an acoustic source embedded in a shear layer using lin- 
earized Euler equations as the far-field equation set are 
encouraging. 

In the present work, the near-field source region so- 
lution including all non-linear flow hydrodynamics is 
obtained through the large-scale equations [8], and is 
matched to the solution of the linearized Euler equations 
governing die acoustic field. The present research com- 
bines the large-scale simulation and linearized Euler 
equations approach into one computer code, resulting in 
the saving of computer CPU time for extension to three 
dimensional acoustic field predictions. 



2. Governing Equations 


The computational domain, figure 1, is split into 
non-linear source generation and linear acoustic propa- 
gation regions, which are governed by the large-scale 
and linearized Euler equations respectively. 


2.1 Large-Scale Equations 


The flow field of a supersonic jet is governed by the 
compressible Navier-Stokes equations, and can be de- 
composed into filtered and residual fields, namely 


f=f+r (1) 

where an overbar denotes the resolved (filtered) field and 
a (“) denotes the unresolved (subgrid) one. The mean of 
the filtered field is the mean of the total field. Upon sub- 
stituting this splitting in the full Navier-Stokes equa- 
tions, the filtered compressible Navier-Stokes equations 
in cylindrical coordinates takes the form 
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Here Q is the unknown vector, F and G are the fluxes in 
the x and r directions, respectively; S is the source term 
that arises in cylindrical polar coordinates; and k is ther- 
mal conductivity. The total enthalpy is I, the total energy 
is E, and Cjj are the viscous stresses. This system of equa- 
tions is coupled with the equation of state for a perfect 


gas. Here, a tilde denotes Favre fitering, 

f = ? (7) 

P 

The unresolved stresses tjj appearing in equations (4-6) 
need to be modeled. 


2£ Linearized Euler Equations 


Starting from the full Navier-Stokes equations in 
conservation form, neglecting viscosity, and linearizing 
about a mean flow (U,V), the axisymmetric linearized 
Euler equations may be written in cylindrical coordinates 
as: 
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and 

p’ = (y-1)[s-(QU+^V) + ip^ + V 2 )] (13) 

In this notation, U and V are the mean axial and radial ve- 
locities. The velocities are normalized by the jet exit cen- 
terline velocity Ue, time by R/U e . density by the mean 
exit centerline value, and pressure by the exit dynamic 
pressure. 


3. Subgrid-Scale Modelling 

The effect of unresolved scales on the resolved ones 
is accounted for through the use of Smagorinsky’s sub- 
grid-scale model [12]. The subgrid-scale turbulence 
stresses are represented as follows]: 
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*lj = -k g (6ij/3) t 2 pv R (s ij -|6 i jS mm ) (14) 

where kg is the kinetic energy of the residual turbulence 
and is neglected with respect to the thermodynamic pres- 
sure. The strain rate of the resolved scale is given by 


ij-alacj axjJ 


(15) 


The summation Smm is zero for incompressible flow, Vr 
is the effective viscosity of the residual field. 
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For the heat equation, Edison [13] proposed the eddy vis- 
cosity model 
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q Ppr t ax k 
where Pr ( is the subgrid-scale turbulent Prandtl number, 
which can be taken as O.S. Smagorinsky’s constant C, in 
equation (16) is 0.1 as was used in previous study [8]. 


4. Numerical Method 


The importance of the dispersion and dissipation of 
a given scheme, used in connection with the computa- 
tional aeroacousdcs, was highlighted by Hardin [14]. 
Both effects are crucial in computational aeroacoustics, 
and can render the computed unsteady part of the solu- 
tion completely unacceptable. As such, high-order accu- 
rate schemes are required for problems in computational 
aeroacoustics. 


A fourth-order accurate in space, second-order accu- 
rate in time scheme is used, which is an extension of the 
McCormack scheme by Gottlieb and Turkel [IS]. This 
scheme has been used extensively by other researchers 
[16-21], to name a few. In this scheme, the operator is 
split into two one-dimensional operators and applied in a 
symmetric way to avoid biasing of the solution. 

5. Boundary Conditions 

Boundary condition is an important issue in the 
computation of jet noise. Proper boundary treatment 
should allow waves to pass through the boundary with- 
out generating reflecting waves. Several boundary treat- 
ments were considered [22], and it was shown that the 
boundary treatments used in [8] were stable, non-reflect- 


ing, and most suitable for the jet computations. The 
present work employs the same boundary treatments as 
in references [8]. The schematic diagram of figure 1 
shows various boundary conditions used at each bound- 
ary. 


In the non-linear acoustic source region, the inflow 
boundary is split into hydrodynamic disturbance and ra- 
diation regimes. A small disturbance, assumed to be 
mainly hydrodynamic in nature, is introduced at the in- 
flow and is specified from the centerline to r/D=2. The 
disturbance is assumed in the form: 

[u* v* p* p] = * (D e‘ (OX ~ } (19) 

Orr-Sommerfeld equation are solved to obtain the com- 
plex wave number a as the eigenvalue corresponding to 
the frequency co and the radial functions d>(r) as the cor- 
responding eigenfunctions. 

O(r) = [h(r),t(r),p(r),j>(r)] (20) 

To obtain the disturbance solution a mean flow must 
be specified. In the present work, the analytical functions 
proposed by Tam and Burton [23] to fit the experimental 
data of Troutt and McLaughlin [24] were used. The mean 
axial velocity is given by 
U = 1 forr<h 

U = exp [-In (2) ^ ) ] f o rr>h (21) 

where b(x) is the half-width of the annular mixing layer, 
and is fitted to the experimental data. The radius of the 
potential core, h(x), is related to b(x) through the conser- 
vation of momentum [23]. 

For the supersonic regime, all characteristics travel 
in the flow direction. Thus the primitive variables are 
given at x=0 as outlined above. In the subsonic regime, 
the following three characteristics are specified accord- 
ing to the linear stability solution: 
p t +pcu t = Cj 

P t -c 2 P = C 2 (22) 

pcv t = C 3 

The fourth characteristic is outgoing and is obtained 
from the interior solution: 

p t -pcu t = C 4 (23) 

The four characteristic equations are then solved togeth- 
er to obtain the time derivatives of the variables, which 
are used to update the solution at the inflow boundary. 

In the radiation regime (r/D>2), the conventional 
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acoustic radiation condition applies: 



where: 

Q = [u,v,p, p] 


R = Jx 2 +t 2 (25) 

r (9 ) = c [|m^-(Im ) 2 

and M is the local Mach number, c is the sonic velocity. 
The spatial derivatives which appear in equation (24) are 
evaluated in an identical manner as the inner flow deriv- 
atives. 

The outflow treatment is based on the asymptotic 
analysis of the linearized equations as given by Tam and 
Webb [25] . The pressure condition is the same as that ob- 
tained by Bayliss and Turkel [26], Enquist and Majda 
[27], and Hariharan and Hagstrom [28], namely: 

p t = -r (0) [-p x + jjP r + j|] (26) 

However, for updating the rest of the primitive variables, 
Tam and Webb have shown that the momentum and con- 
tinuity equations should be used to account for the pres- 
ence of entropy and vorticity waves at the outflow 
boundary. 

For the outflow regime of large radius with the local 
Mach number less than 0.01, the outflow condition is re- 
placed by the conventional acoustic radiation condition, 
equation (24). 

The results presented herein are for axisymmetric 
disturbance, for which the boundary condition at r=0 can 
be stated as 

s[«.p.p] = °^ = ° <*> 

The centerline treatment for non-axisymmetric distur- 
bances is not obvious, and was addressed by Shih et al 
CT. 

In the linear acoustic region, radiation boundary 
condition of equation (24) is used since only acoustic 
wave is of significance in this region. 

6. Matching 

In the present work, this scheme is applied to both 
the large-scale equations and linearized Euler equations. 
A schematic diagram of the computational domain is 


presented in Figure 1. No boundary condition is needed 
in the matching region. The scheme requires fluxes at 
two grid points outside of the domain to update variables 
at the boundary. An overlapping region, consisting of 4 
grid points in the radial direction, is made to allow vari- 
ables at interior points to be passed between two do- 
mains. The variables of the two solution domains are the 
filtered quantities and the fluctuation quantities. The 
mean of the filtered field is the mean of the total field. 
Therefore, it is straightforward to set up a matching 
scheme at die interface. At every time step, all linear 
acoustic region variables in the overlap region are con- 
verted to non-linear region variables and used to evaluate 
the fluxes in the non-linear region. Likewise at every 
time step, all non-linear region variables are converted to 
linear acoustic region variables and used to evaluate the 
fluxes in the acoustic region. Freund et al. [10] pointed 
out that grid stretching and interpolation induce error in 
the calculation, though reducing the computing expense. 
The present simulation uses identical uniform grids in 
the overlapping region to avoid such generation of spuri- 
ous waves. 

7. Results and Discussion 

The numerical simulation was conducted for a cold, 
nearly perfectly expanded axisymmetric supersonic jet 
of Mach number 2.1. The test case has been simulated 
previously [8] using large-scale equations. The total tem- 
perature of the jet is 294° K, and the jet exit pressure is 
0.0515 atm. The Reynolds number based on exit condi- 
tions is approximately 70000. In the present calculation, 
the jet is excited at a Strouhal number of 0.2 with the 
Strouhal number defined as St=fD/U e , where D is the 
nozzle exit diameter and U c is the jet exit centerline ve- 
locity. 

The computational domain for this problem extends 
axially from x/D=2.5 to x/D=35, and radially from cen- 
terline to r/D=16, as shown in figure 1. Due to the steep 
mean flow gradient encountered at the jet exit, the com- 
putational grid was begun at an axial distance x/D=23 
from die actual jet exit The computational grid consists 
of 391 equally spaced points in the axial direction. In the 
radial direction, 150 points are used and stretched be- 
tween centerline and r/D=2.5 with concentration of grid 
points around r/D=05. Between t/D=2S and r/D=16, 
130 equally spaced points are used with a spacing equal 
to that of the last stretched points. The interface between 
the non-linear and linear acoustic regions was chosen at 
r/D=5, above which the mean flow is uniform throughout 
the domain. 
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Figure 2 shows the instantaneous distribution of p, 
p.u and vat t=150, where t is the characteristic time de- 
fined as the ratio of nozzle exit radius to the jet exit cen- 
terline velocity. The horizontal line at t/D=5 in the figure 
Indicates the boundary between the two domains. One 
can see smooth solutions across the ma tc h i n g region, al- 
though some kinks are observed in the axial velocity dis- 
tribution. The disturbances generated in the shear layer 
propagate through the matching region to the far field 
without any distortion. The wave-like nature of the flow 
field is evident and the solution is clean from boundary 
reflections. 

Figure 3 shows contours of the root-mean-square 
values of the pressure and axial momentum distributions. 
The contour levels are 0.00035 and 0.0025 for die root- 
mean-square values of pressure and axial momentum re- 
spectively. No spurious disturbances are generated in the 
matching region, and the preferred forward emission is 
clearly shown in the figure. 

Figures 4 shows the axial development of the axial 
momentum fluctuations on the nozzle lip line. The com- 
puted results obtained using the large-scale equations 
throughout the computational domain [8] are also shown 
in this figure. The present zonal approach gives die de- 
velopment of axial momentum fluctuations identical to 
those obtained in reference [8], as was expected. Figure 
5 shows the axial development of the root-mean-square 
values of the pressure on the nozzle lip line. Again, the 
result of zonal approach is identical to those in reference 
18]. 

Figure 6 compares the amplitude of pressure distur- 
bance at r/D=12 along the axial direction at St=0.2 with 
the results from reference [8]. Due to the use of linear- 
ized Euler equations in this region (no viscous dissipa- 
tion), slightly higher peak was obtained by die present 
approach. The overall distribution is consistent with the 
previous large-scale computed results. 

Figure 7 presents the pressure spectra at x/D=25, r/ 
D=12 for both zonal and large-scale approach. The am- 
plitude at fundamental frequency, StsO.2, is slightly 
higher for the zonal approach, while becomes lower at 
successive harmonics, than the large-scale simulation re- 
sults. Previous linearized Euler computation along [29] 
did not show any generation of successive harmonics ex- 
cept for the fundamental one. 

Figure 8 shows the sound pressure level distribution 
in the far field for the present and previous calculations. 
It is seen that the computed patterns of the sound pres- 
sure level contours are consistent with the previous work 


[8]. As pointed out in reference [8], the computed results 
show a downstream shift of the lobes when compared to 
the experimental measurements. This can be explained 
that the calculations are for axisymmetric case, while the 
measurements showed that excited motion of die jet is 
dominated by die first helical mode in addition to the ax- 
isymmetric mode [24]. 

Figure 9 shows the calculated sound field directivity 
at a circle of radius 24D with center at the jet exit center- 
line. The angle is measured from the jet exit centerline. 
The calculated peak occurs around 15 degrees, which is 
the same as the one predicted by previous work [8]. 

The solution sensitivity to matching boundary loca- 
tion, which was chosen at r/D=5 and 7.5 respectively, is 
also investigated. Figure 10 shows die comparison of the 
axial distribution of root mean square values of pressure 
at r/D= 5 and 7.5. As one can see that the solution is in- 
dependent of the matching boundary location, which is 
denoted by RM in the figure. 

8. Conclusions 

A zonal approach for direct computation of sound 
generation and propagation from a supersonic jet is in- 
vestigated. The computational domain is split into a non- 
linear acoustic source generation regime and a linear 
acoustic wave propagation regime. The unsteady flow in 
the non-linear acoustic source region is governed by the 
large-scale equations, which are the filtered compress- 
ible Navier-Stokes equations. The linearized Euler equa- 
tions are used to describe the sound wave propagation in 
the linear acoustic region. The computed results show 
that no spurious wave is generated in the matching re- 
gion and the computational cost is reduced by 30% when 
compared with the direct simulation using large-scale 
equations alone. This substantial reduction in computa- 
tional cost becomes significant when the approach is ex- 
tended to three dimensional simulations. 
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Fig. 1 Computational Domain 
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Fig. 2 Instantaneous distributions of (a) density (b) axial velocity (c) pressure (d) radial velocity, 
at t=150. 



Fig. 3 Contours of root mean square values of (a) pressure (b) axial momentum. 
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Fig. 4 Root mean square values of axial 
momentum at r/D=0.5. 


Fig. 5 Root mean square value of pressure at 
r/D=0.5. 


Strouhal number 


Fig. 6 Amplitude of pressure fluctuations at r/ Fig. 7 Pressure spectra at x/D=25, r/D=12. 


Fig. 8a Sound pressure level contours for 
large-scale calculations. 


Fig. 8b Sound pressure level contours for laige 
scale and linearized Euler calculations. 








Fig. 9 Directivity of jet noise at R/D=24. 



Fig. 10a Comparison of root mean square value 
of pressure at r/D=5. 



Fig. 10b Comparison of root mean square value 
of pressure at r/D=7.5. 
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